Introgression between <i>Betula tianshanica</i> and <i>Betula microphylla</i> and its implications for conservation

In rapidly changing environments species conservation can be hindered by uncertainties in distinguishing closely related species. Cryptic ongoing hybridization add further uncertainty and could beneficial or destructive. Here, we show that a declining birch tree is hybridizing with more widespread relative the Junggar basin, NW China, their hybrids have been previously named as rare sub-species. Given numbers of this tree, suggest effort should aim to slow rate anthropogenic habitat loss at hybrid zone preserve its pure populations away from zone. Hybridization play an important role plant evolution introducing novel alleles (Goulet et al., 2017) leading speciation (Abbott 2013; Whitney 2010). However, also lead significant challenges where it involves common (Balao 2015), native invasive (Bleeker 2007), domesticated wild ancestor (Wayne & Shaffer, 2016). potentially reduce fitness via genetic demographic swamping, which may introgress maladaptive (Rhymer Simberloff, 1996; Todesco Such situation worse for wind-pollinated because pollen transfer long distances (Hjelmroos, 1991; Skjøth far beyond current species’ distribution. The genus Betula excellent model assess consequences among consists approximately 65 subspecies (Ashburner McAllister, 2016; Wang 2016), extensive introgression documented based on morphology, cytogenetics, molecular markers (Anamthawat-Jónsson Thórsson, 2003; Anamthawat-Jónsson Tómasson, 1990; Bona 2018; Eidesen 2015; Hu 2019; Tsuda 2017). are pioneers, colonizing new habitats providing shelter other trees. fact often co-occur means frequent them Some restricted distributions, such B. nana United Kingdom, possibly due historical fragmentation (Borrell 2014; Zohren Several others endangered, megrelica Georgia calcicola southwest China study, investigate possible decline regional endemic microphylla, verge extinction rapid climate change human activities (Zhang Kong, 2019) (Figure 1). microphylla generally grows wetlands within arid semi-arid areas evaporation much higher than precipitation. Many seem gone extinct during past decades basin. For example, Caotanhu wetland, was historically abundant existed 2011, but population had disappeared 2018 (pers. comm.). A decrease level groundwater activities, resulting high salinity, thought responsible wetland (Li Zhang 2019). has overlapping range tianshanica both tetraploids mainly distributed Altay Mountains (to north Xinjiang province). Historical occurrence data used occur Tianshan common. two tetraploid very likely sister species, forming morphological continuum no clear boundaries S1), hindering efforts. province, several varieties proposed minor differences, var. ebinurica, according locality discovered (Yang 2006). studies suggested these need conserved urgently (Duan Yang systematic classification available guide conservation. ranges overlap platyphylla latter diploid, any unlikely viable. SSRs restriction site-associated DNA sequencing (RAD-seq) determine whether microphylla: ebinurica represent distinctive cluster deserve taxonomic recognition consideration addition, inter-specific admixture occurs between microphylla. We characterized local diversity check our focal Our specific questions follows: (a) Do tianshanica, form distinct clusters? (b) Does platyphylla? (c) If occurs, does pose threat microphylla? To end, genotyped 265 individuals collected 10 15 microsatellite conducted RAD-seq 49 four five randomly chosen each population. 42 obtained three using markers. Populations were located records samples July August 2018. All appeared natural. Within population, separated least 20 meters avoid sampling same clone. Sampling locations recorded GPS system (UniStrong). herbarium specimen created individual dried press. total Mountains, Mountains. Based morphology location, SLM, EYQ, EYQE treated HBK, TC, JMN, YBW, XY ABH MGH basin tentatively defined 2006) comm.), respectively. 23, 17, XEXL, XY, total, 83 53 129 potential distribution platyphylla, sourced coordinate points fieldwork literature. 23 26 being 5 km apart other. Nineteen bioclimatic layers downloaded WorldClim database (http://www.worldclim.org/) 2.5 arc-min resolution, LGM (CCSM v4 scenario), MID (CNRM-CM5 (the period 1960–1990), future periods (2050 2070, MIROC-ESM scenario). Highly correlated variables (correlation coefficient > 0.7) identified “pairs” function R package “raster” version 3.1–5 (Hijmans, removed overfitting, six retained analysis. These low correlation BIO1 (annual mean temperature), BIO2 (mean diurnal range), BIO4 (temperature seasonality), BIO9 temperature driest quarter), BIO12 precipitation), BIO15 (precipitation seasonality). An ecological niche (ENM) developed MaxEnt (Phillips performed 50 subsampled replicate runs 25% observations cross-validation. Models evaluated receiver operating characteristic (ROC) area under curve (AUC). AUC value greater 0.9 indicates good prediction, while 0.7 predicted better random (Swets, 1988). compute calculated Schoener's D Hellinger distance I (Schoener, 1970; Warren 2008) “raster.overlap” “ENMTools” 0.2 (https://github.com/danlwarren/ENMTools). Five intact mature leaves pressed selected morphometric scanned 1,325 Hewlett-Packard printer (LaserJet Pro MFP M128fn) resolution 600 d.p.i. 13 landmarks leaf described previous study (Hu Liu 2018). Briefly, coordinates thirteen digitized converted into configuration cartesian pairs (x, y) program ImageJ (Abràmoff 2004). Principal component analyses (PCAs) normalized matrix specimen, respectively, MORPHOJ (Klingenberg, 2011). Genomic isolated cambium tissue following modified 2 × CTAB (cetyltrimethylammonium bromide) protocol (Wang 2013). quality genomic assessed 1.0% agarose gels then diluted concentration 10–20 ng/µl genotyping. Fifteen loci pendula (Kulju 2004), pubescens ssp. tortuosa (Truong 2005), maximowicziana (Tsuda 2009), japonica (Wu 2002) genotyping (Table S1). 5′ terminus forward primers labeled FAM, HEX, TAM fluorescent probes. Each locus amplified individually prior artificially combined multiplexes order errors caused size overlapping, length differences dye. PCR provided supplementary data. Microsatellite scored software GENEMARKER 2.4.0 (Softgenetics) checked manually. Individuals missing excluded, 307 final dataset. delivered KEGENE Company library preparation enzymes PstI MspI digest DNA. Libraries sequenced Illumina HiSeq 2500 pair-end 150-bp sequencing. Adapter sequences custom Perl script. raw cleaned process_radtags module Stacks v2.2 (Catchen 2013) paired-end mode. sliding-window step window 15% read. Reads required below 30 unpaired reads discarded. Clean aligned whole genome (Salojärvi BWA-MEM v.0.7.17-r1188 algorithm BWA (v0.7.17) default parameters Durbin, 2009). non-specific matches One excluded subsequent analyses, samples. Alignments sequence alignment map (SAM) format sorted, indexed binary (BAM) files (SAMtools v1.8) MarkDuplicates tool Genome Analysis Tool Kit (GATK) (v 4.1.4.) mark duplicates (DePristo 2011; McKenna HaplotypeCaller GATK call genotypes sample (v4.1.4), generating intermediate gVCF files. CombineGVCFs merge VCF file, input file joint GenotypeGVCFs tool. SNPs filtered mapping (MQ) threshold 20, variant confidence (QUAL) 30, QUAL score (QD) 2, maximum symmetric odds ratio (SOR) 3, minimum depth (DP) 200, probability strand bias (FS) excess heterozygosity (ExcessHet) 60 54.69, Z-score read qualities (MQRankSum) position (ReadPosRankSum) −12.5 −8.0, SelectVariants filtering applied select present 50% BCFtools v1.8 remove kb r2 .5 linkage disequilibrium allele frequency (MAF) <0.01 (Li, deposited NCBI-Sequence Read Archive (SRA) repository BioProjectID PRJNA679411. (PCO) analysis SSR POLYSAT (Clark Jasieniuk, 2011) implemented 3.6.1 (R Core Team, pairwise Bruvo al. (2004). nucleotide SNPs, principal (PCA) carried out “adegenet” 2.1.1 (Jombart, 2008), setting ploidy four. Both datasets analyzed STRUCTURE 2.3.4 (Pritchard 2000) identify most number clusters (K) Ten replicates computational limitations. 1,000,000 iterations burn-in 100,000 run K 1 generations (five replicates) SNPs. model, assumption frequencies populations. assigned highest membership averaged over ten independent runs. make biological sense estimated “Evanno test” (Evanno 2005) Structure Harvester (Earl vonHoldt, 2012) “Thermo-dynamic Integration” (TI) method (Verity Nichols, 2016) As did not use TI designed thousands Replicate grouped symmetrical similarity >0.9 Greedy CLUMPP (Jakobsson Rosenberg, 2007) visualized DISTRUCT 1.1 (Rosenberg, Gene (Nei, 1987) allelic richness (El Mousadik Petit, 1996) FSTAT 2.9.4 (Goudet, 1995) diploid (2017) (2019). plotted derived against optimal value. Then, Mann–Whitney U-test test difference different marker sets “wilcox. test.” statistical if latitude performing mixed effects “lme” (Pinheiro Bates, constructed models well ranged 0.932 0.940 SD 0.025 0.029 timescales satisfactory 0.813 0.839 0.044 0.065. environmental predictors annual precipitation seasonality Suitable limited exist western part outside Xinjiang, whereas suitable west 2). During Holocene, along environs. seems stable Holocene larger 2050 shift northwards. contracted since become fragmented 2070 conservative (0.40) showed little (I = 0.036, 0.001) 0.229, 0.115) 0.255, 0.216). Range 0.362, 0.303) 0.349, 0.298) 2070. consistently broader all PCA revealed single cluster. PC1 PC2 explained 30.3% 18.7% variation leaves, respectively 3), summarized 40.9% 23.4% trees, S2). PCO Bruvo's 4) accounted 26.7% variation. differentiated formed separate closer together, some bridging gap them. Consistent results, clusters. “TI” 3 value, ΔK criterion failed S3a b). Cluster I, shaded orange Figure 5, included EYQE, ABH, MGH, correspond tianshanica. II, blue III, green corresponded 5a). relatively there seemed one early-generation SLM. TC HBK YBW admixture. Limited detected sympatric found regarded containing sub-species result suggests level, marginally per 6,436,408 92,924,972 trimmed 5,827,435 84,297,044 5,165,913 78,364,760 mapped reference sample. mappings resulted 64.1% 93.7% individual. Over 90% successfully After (> million reads, only SNVs, coverage <50% data), 386,363 sites 661 individuals. basis variants individuals, 16,110 biallelic, 493 triallelic, tetra-allelic. S4) genotype calls 16,616 indicated corresponding 5b) supported S5). EYQ levels 0.04% 0.2%, S6). More values 8.7%, 5.1%, 12.1%, 5b, S7). 5b c). 23.9% 18.1% lower northerly 14.0%, 1.8%, 1.8% c, S7a). At K, subdivision occurred populations, JMN S7b). distinguished S7a (p .1386, U test) yielded similar estimates significantly 6a, p .0002, test). Similarly, SLM after excluding SLM020 6b, .004, Levels decreased south .127, S8a) .057, S8b). strong evidence they genetically entities despite referred history plausible geographically main Previous research noted resembled On findings, sub-species, SNP apparent vice versa. tended 5a b), pattern comparison Bradbury (2015) (Bohling individual, namely, SLM020, 51.8% b; This smaller SSRs, less representative genomes. Introgression (SLM MGH) reported past. proximity without Thus, center hybridization. presence tianshanica-derived involving unexpectedly flow undiscovered region. It selection adaptive introgression, adapting Another explanation incomplete lineage sorting (ILS) geographic signals (Barton, 2001) ILS making contribution patterns face threats loss. ENM results predict northward expansion future, expect even currently. field observations, usually growing riversides mountain slopes adjacent rivers so advance helped seed dispersal flowing observed seedlings samples, including provide degradation. Without recruitment cannot survive. northern regions grow wetlands, surrounded vastly inhospitable lands, thus limiting successful seeds. Overgrazing crop plantation modifying greatly, posing unprecedented risk existence Shihezi city facing activities. implications management viewed ways, depending upon seen creative destructive process system. will possibilities adaptation Vallejo-Marín Hiscock, ensuring persistence degeneration habitats. view leads process-based approach conserving evolutionary processes generate biodiversity (Ennos 2012). taxonomical complex taxa, Sorbus (Hamston Robertson Under view, desirable conserve parent hybrids, those known (cf. Duan Li hand, argued system, contributing greatly success simply generation invades northwards, rather superior climates. case, desire measures include removal now, course would wise seeks enhance eliminate Clearly, needed understand gene detrimental context change. beneficial, persist expand. immediate stop maintain perhaps propagating planting becoming thank Dr. Laura Kelly Royal Botanic Gardens Kew valuable comments manuscript Jihong Huang Research Institute Forest Ecology, Environment Protection, Chinese Academy Forestry sharing us. work funded National Natural Science Foundation (31770230 31600295) Funds Shandong ‘Double Tops' Program (SYL2017XTTD13). DH conceived project. DH, NW, LW JD, LW, ZL, FW laboratory work. JD results. wrote initial draft. RJAB JSB substantial thereafter. Please note: publisher content functionality supporting information supplied authors. Any queries (other content) directed author article.